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Emission at far-infrared wavelengths makes up a significant fraction of the total light detected from galaxies over the age of Universe. Herschel 
provides an opportunity for studying galaxies at the peak wavelength of their emission. Our aim is to provide a benchmark for models of galaxy 
population evolution and to test pre-existing models of galaxies. With the Herschel Multi-tiered Extra-galactic survey, HerMES, we have observed 
a number of fields of different areas and sensitivity using the SPIRE instrument on Herschel. We have determined the number counts of galaxies 
down to ~ 20 mJy. Our constraints from directly counting galaxies are consistent with, though more precise than, estimates from the BLAST 
fluctuation analysis. We have found a steep rise in the Euclidean normalised counts < 100 mJy. We have directly resolved ~ 15% of the infrared 
extra-galactic background at the wavelength near where it peaks. 

Key words, galaxies: evolution, submillimeter: galaxies, surveys 



1. Introduction 

The statistical properties of galaxy populations are important 
probes for understanding the evolution of galaxies. The most 
basic statistic of galaxy populations is the number counts i.e. 
the number density of galaxies as a function of flux. The first 
strong evidence for cosmological evolution came through study- 
ing number counts of radio galaxies (e.g. Longair 1966). 

The number counts at far-infrared and sub-mm wavelengths 
are well known to exhibit strong evolution, e.g. from IRAS 
(jOliver et al.||1992| and references therein), ISO ( jOliver et al.| 



BLAST, that count models have been probed using fluctuation 
techniques (Patanchon et al. 2009) or directly (Bethermin et al. 
2010) at the wavelength where the background peaks. 

Far-infrared and sub-mm counts and background measure- 
ments have been modelled phenomenologically with strongly 
evolving populations (see Section [4}. Physical models (e.g. so- 
called semi-analytic models) struggle to explain these counts 
and solutions include altering the initial mass function (e.g. 
Bau gh et al.||2005[) or exploiting AGN/Supernovae feedback 



(e.g. |Granato et al.|2004) . 
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> |2002| jH eraude au et al.||2004| and ref erences therein) Spitzer a primary goal of Herschel (|Pilbratt et al7||20TQ> is to ex 



^-j ( |Shupe et al.||2008[ |Frayer et al.||2009| and references therein) 



and ground-based sub-mm surveys (Maloney et al. 2005 ; Coppin 
et al.|2006HKhan et al.|2007i|Greve et al.|2008HWeifi et al.|20091 



Scott et al.|2010| and references therein). 

These results are underlined by the discovery of a signifi- 



cant extragalactic infrared background (Puget et al. 1996 Fixsen 
|et al.| 1998; Lagach e~et al.|1999] l. The background measures the 
flux weighted integral of the number counts over all redshifts 
plus any diffuse cosmological component. This indicates that as 
much energy is received from galaxies after being reprocessed 
through dust as is received directly. It is only very recently, using 



plore the evolution of obscured galaxies. Herschel opens up a 
huge region of new parameter space of surveys in area, depth 
and wavelength. 

The Herschel multi-tiered extragalactic survey (HerMESQ 
Oliver et al. 2Q10|> is the largest project being undertaken by 



Herschel and consists of a survey of many well-studied extra- 
galactic fields (totalling ~ 70deg 2 ) at various depths. This let- 
ter is the first number count analysis from the HerMES Science 
Demonsnation Phase SPIRE data. Even these preliminary re- 
sults will be able to eliminate some existing models and provide 
a benchmark on which future models can be tested. 



* Herschel is an ESA space observatory with science instruments 
provided by European-led Principal Investigator consortia and with im- 
portant participation from NASA. 



hermes.sussex.ac.uk 
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2. SPIRE Data 

2.1. Science Demonstration Phase Observations 

The observations described here were carried out on the 
Herschel Space Observatory (Pilbratt et al. 2010) using the 
Spectral and Photometric Imaging Receiver (SPIRE). The 
SPIRE instrument, its in-orbit performance, and its scientific ca- 
pabilities are described by Griffin et al. (2010), and the SPIRE 
astronomical calibration methods and accuracy are outlined in 



Swinyard et al.| ( |2010] >. They were undertaken as part of the 
HerMES programe during the Science Demonstration Phase be- 
tween 12-Sep-2009 and 25-Oct-2009 under the proposal identi- 
fication SDP_soliver_3. The fields and observations are sum- 
marised in Table Q] 



2.2. SPIRE Catalogue Data Processing 

For this paper it is sufficient to note that the same source de- 
tection method is applied to the simulations as to the real data 
so we sketch the details only briefly. The single-band SPIRE 
catalogues have been extracted from the maps using a version 
of the sussExtractor method (Savage & Oliver 2007) as imple- 
mented in hipe (Ott 2010). The processing of the SPIRE data is 
summarised here, with details of the approach given by Smith| 
et al. ( 2010| l. Calibrated timelines were created using HIPE de- 
velopment version 2.0.905, with a fix applied to the astrometry 
(included in more recent versions of the pipeline), with newer 
calibration files (beam-steering mirror calibration version 2, flux 
conversion version 2.3 and temperature drift correction version 
2.3.2) and with a median and linear slope subtracted from each 
timeline. The default hipe naive map-maker was then used to 
create maps, which were given a zero mean. The shallow fields 
(Lockman-S WIRE and FLS) were smoothed with a point-source 
optimised filter (see Smith et al. 2010 for details). Peaks in the 
map were identified and the flux was estimated based on an as- 
sumed (Gaussian) profile for a point source, through a weighted 
sum of the map pixels close to the centre of the source. This 
filtering means we underestimate the flux of extended sources. 
The SPIRE Catalogue (SCAT) processing is assessed by inject- 
ing synthetic sources on a grid into the real maps. We then run 
the SCAT source extraction pipeline on these maps and claim 
success if the closest detection to the injected source is within a 
search radius of fwhm of the beam and has a flux within a factor 



of two of the injected flux (see Smith et al. 2010 for more de- 
tails). The resulting 50% completeness estimated in this way is 
tabulated in TableQ]but is not used to assess the counts. 



3. Method 

The Herschel beam is broad compared with the number density 



of sources, i.e. the maps are confused. Nguyen etal. ( 2010) mea- 
sure a variance in nominal map pixels due to confused sources 
finding cr con f 5.8, 6.3 and 6.8 mJy/beam. This confusion means 
care has to be taken in the estimation of number counts. Our 
technique follows the standard approach for sub-mm surveys, 
correcting for flux boosting and incompleteness. 

We determined the false detection rate by applying the 
source extraction on maps obtained from the difference between 
two independent observations of the same field. These maps are 
expected to have zero mean, no sources, but similar noise prop- 
erties to mean maps. We thus estimate that the reliability for the 
samples in this paper is better than 97% for all fields and bands. 



A source we measure to have flux, S m , and noise, cr m , is 
more likely to be a dimmer source on top of a positive noise 
fluctuation than the converse; this is known as flux boosting. 
We follow the Bayesian method of Crawf ord et al.| ( |2009| l for 
estimating fluxes of individual sources ("de-boosting"). We es- 
timate the posterior probability distribution of the true flux of 
the source (Si) that contributed the most flux to a given detec- 
tion. Note that it is similar to the now-standard flux de-boosting 
method (e.g. Copp in et al.||2006] > but with an additional expo- 
nential suppression term at low intrinsic flux. We derive counts 
by randomly sampling the posterior distribution ten thousand 
times. The flux de-boosting procedure has some dependency on 
the choice of prior number count model and so these samples 
are drawn from distributions produced for the full range of mod- 
els discussed in Section HI These samples provide a direct esti- 
mate of the confidence region of our counts. For faint sources, 
the posterior probability function rises beyond the sampled flux 
range (S m /5 < S, < 5S m ) so that the deboosted flux is highly 
uncertain. In those cases we flag the deboosted flux as "bad" 
and the derived number counts at the flux level where deboosted 
fluxes are flagged "bad' are unreliable. We exclude count bins 
in which the fraction of "bad" sources is > 20%. We also es- 
timate the uncertainty in this by looking at the variation in de- 
rived counts from the range of models. Errors are included in the 
plots. The flux de-boosting procedure assumes no clustering. 
Clustering will affect this and will be addressed in a later paper. 

We estimate the incompleteness in the whole process by 
running full simulations. We have constructed input maps var- 
ious from various number count models (Pearson & Khan 2009 



Lagache et al. 2004 ; Patanchon et al. 2009 and models from Xu 
et al. (2003) and Lacey et al. (2010)). These input maps are then 



9; 
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processed by the SPIRE Photometer Simulator (sps, Sibthorpe 
et al. 2009) for observational programmes exactly the same as 
the real data. The timeline output of the sps, map-making and 
source extraction are then processed in the same way as the real 
data, including the flux-deboosting. The ratio of input to output 
counts gives us the completeness with the standard deviation be- 
tween models providing an estimate of an error in that estimate. 

4. Results 

The results are presented in Table [2] and in Figure [T] with 
Euclidian normalization. There are several sources of uncertain- 
ties in the number counts: Poisson noise from the raw counts; 
"sampling variance" due to additional fluctuations from real 
large-scale structure; additional Poisson noise from the sampling 
of the posterior flux distribution and systematic errors from the 
corrections and assumptions about priors and the effect of clus- 
tering on the de-boosting. We measure the standard deviation of 
counts between fields which includes Poisson errors and some 
of the other systematic errors. The errors plotted are the field 
to field variations (or the Poisson errors if larger) with the er- 
rors from flux-boosting and completeness corrections added in 
quadrature. 

We see approximately flat counts for S > 100 mJy and then 
a steep rise. There is flattening to about 20 mJy. We find very 
good agreement with the number counts estimated from a P (D) 
fluctuation analysis of the BLAST maps (Patanc hon et al.|2 009 ). 

We have also estimated, but do not show, the integral counts. 
The flux density at which the integral source counts reach 1 
source per 40 beams (with beams defined as 3.87 x 10~ 5 , 7.28 x 
10- 5 ,1.48xl0~ 4 deg 2 )is 18.7 + 1.2, 18.4+1.1 and 13.2±1.0mJy 
at 250, 350 and 500 fim respectively ( N.B. these fluxes are 
slightly below our secure estimation of counts) . Likewise the 



S. J. Oliver et al.: HerMES: SPIRE galaxy number counts at 250, 350 and 500 yum 



Table 1. HerMES SPIRE SDP Observations. PACS observations will be discussed in Aussel et al. (2010 1. Size is approximate 
extent of region with typical coverage. Roll angle is measured East of North. Repeats is total number of pairs of scans in both A and 
B directions. ?aor is total time in execution of the observations. (N ssm p) is the mean number of bolometer samples per pixel in the 
same typical-coverage region of the 250 pm map (6" x 6"pixels). S 50 % m is the flux density at which 50% of sources injected into 
the 250 pm map are faithfully recovered. 



Name 


Size 


RA 


Dec 


Roll 


Mode 


Scan 


Repeats 


?AOR 


V'samp/ 


„ 250/jm 
^ 50% 


r, 35U/jm 
^50% 


„ 500pm 
^50% 






r 


r 


r 




Rate 




fttt 




/mJy 


/mJy 


/mJy 


A2218 


9'x9' 


248.98 


66.22 


217 


Lrg. Map 


30"/s 


100 


9.2 


1622 


13.8 


16.0 


15.1 


FLS 


155' x 135' 


258.97 


59.39 


185 


Parallel 


20"/s 


1 


16.8 


30 


17.5 


18.9 


21.4 


Lockman-North 


35' x 35' 


161.50 


59.02 


91 


Lrg. Map 


30"/s 


7 


3.9 


117 


13.7 


16.5 


16.0 


Lockman-SWIRE 


218' x 218' 


162.00 


58.11 


92 


Lrg. Map 


60"/s 


2 


13.4 


16 


25.7 


27.5 


33.4 


GOODS-N 


30' x 30' 


189.23 


62.24 


132 


Lrg. Map 


30"/s 


30 


13.5 


501 


12.0 


13.7 


12.8 



number density at 100 mJy is 12.8 + 3.5, 3.7 ± 0.4 and 0.8 ± 
0. 1 deg~ . These last measurements alone will be sufficient to 
rule out many models. 

Since the first ISO results, many empirical models have been 
developed to predict and interpret the numbers and luminosi- 
ties of IR galaxies as a function of redshift. Empirical models 
are based on a similar philosophy. The spectral energy distri- 
butions of different galaxy populations are fixed and the mid- 
IR, far-IR and submm data are used to constrain the luminosity 
function evolution. Current limits come from the mid-IR, far- 
IR and submm number counts, redshift distributions, luminosity 
functions, and cosmic IR background. Models all agree on the 
general trends, with a very strong evolution of the bright-end 
(> 10 11 L Q ) of the luminosity function and they yield approx- 
imately the same comoving number density of infrared lumi- 
nous galaxies as a function of redshift. We compare the num- 
ber counts with eight models, one pre-Spitzer ( |Xu et al. 2003), 
two based on the ISO, SCUBA and Spitzer first results (Lagache 
et al.|2004 Negrello et al.|2007 ) and 5 being more constrained 



by de ep Spitzer, SCUBA, AzTEC, and recent BLAST observa- 
tions (|Le Borgne et al.||2009"l [Pearson & Khan||2009"l |Rowan-| 
|Robin son 2009; Va liante et al . 2009, Franceschini et al.||2010) 
The differences between the models are in several details, dif- 
ferent assumptions leading sometimes to equally good fits to the 
current data. For example, Valiante et al. ( 2009) conclude that it 
is necessary to introduce both an evolution in the AGN contri- 
bution and an evolution in the luminosity-temperature relation, 
while Franceschini et al. (2010! reproduce the current data with 
only 4 galaxy populations and only one template for each popu- 
lation. We also compare with two semi-analytic models those of 



Lacey et al.| ( (2010| > and | Wilman etaT[ ( [2010) 



Comparison with SPIRE number counts shows that many 
models cannot fit the bright end (> 100 mJy). Exceptions are 
the models of |Negrello et al.| p007), |Valiante et al.| ( J2009) , 
|Franceschini et al.|([2010[> a nd |Pearson & Kh an (2009]TOf these 
only Valiante et al.| (|2009[) c an fit the rise from (20 < S < 
100) mJy. The |Valiante et al.| ( |2009| l model has "cooler" spectral 
energy distributions at higher redshift. However, increasing the 
number of higher redshift galaxies would have a similar effect 
on the counts so it would be premature to assume the spectral 
energy distributions need revision. 

We have also calculated the contribution of the resolved 
sources to the background intensity as a function of flux (shown 
in Figure 0). At the (40 beams) -1 depth we resolve 1.73 ± 
0.33, 0.63 + 0.18, 0.15 ±0.07 nWm 2 or 15, 10, 6% of the nom- 



inal measured values at 250, 350 and 500 pm (Lagache et al. 
[19991 . 

Future work will provide more detailed constraints at the 
fainter limits. This will include a P(D) analysis ( |Glenn et al.| 



Table 2. Table of HerMES SDP number counts at 250, 350 and 
500pm. Bin limits and Euclidian weighted central fluxes are 
given in mJy. Counts are in /sr~' Jy~' . Errors 1 - 4 are fractional 
errors in percentages arising from: flux-deboosting, complete- 
ness corrections, Poisson statistics and field-field variations re- 
spectively. 





Binl 


Bin 2 


Bin 3 


Bin 4 


Bin 5 


Bin 6 


^ min 


20 


29 


51 


69 


111 


289 


•J max 


29 


51 


69 


111 


289 


511 


u euc 


23.8 


37.5 


58.9 


85.9 


166.2 


374.1 








250yum 








dN 


2.0xl0 8 


6.4xl0 7 


1.2x1 7 


3.1xl0 6 


2.1xl0 5 


1.7xl0 4 


erri 


38 


22 


28 


35 


23 


6 


erri 


10 


6 


7 


6 


4 


19 


err3 


4 


3 


7 


4 


7 


23 


err 4 


11 


19 


30 


28 


30 


8 








350yum 








dN 


1.1x10 s 


3.5xl0 7 


5.3xl0 6 


l.lxlO 6 


6.2xl0 4 


4.7xl0 3 


erri 


49 


34 


44 


56 


25 


129 


erri 


18 


14 


17 


6 


5 


12 


err 3 


7 


4 


10 


7 


15 


43 


err 4 


23 


13 


33 


8 


12 


64 








500yum 








dN 


3.6xl0 7 


l.lxlO 7 


1.6xl0 6 


2.3xl0 5 


1.3xl0 4 


1.3xl0 3 


cri'i 


83 


50 


62 


56 


45 





erro 


31 


18 


25 


15 


18 


7 


eiT3 


10 


6 


18 


14 


33 


50 


err 4 


5 


17 


48 


27 


20 






2010) and counts from catalogues extracted at known Spitzer 
source positions. 

5. Conclusions 

We present the first SPIRE number count analysis of resolved 
sources, conservatively within the limit of Herschel confusion. 
We have measured counts which resolve around 15% of the 
infrared background at 250 /mi. We see a very steep rise in 
the counts from 100 to 20 mJy in 250, 350 and 500 pm. Few 
models have quite such a steep rise. Many models fail at the 
bright counts > 100 mJy. This may suggest that models need a 
wider variety or evolution of the spectral energy distributions or 
changes in the redshift distributions. Future work is required to 
accurately constrain the fainter ends < 20 mJy where confusion 
is a serious challenge. 
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Fig. 1. Number counts obtained from HerMES source cata- 
logues. Filled circles are the mean number counts averaged over 
the following fields. GOODS-N & Lockman-North (faintest five 
bins only) and FLS & Lockman-S WIRE (brightest six bins only) 
with flux-deboosting, completeness corrections and field-field 
error bars. Model fit to fluctuations of BLAST maps (omitting 
upper-limits, Patanc hon et al.|2009|) — shaded region; BLAST 
resolved counts (Bethermin et al.|2010J> — open triangles; Khan 



et al. ( 2007 1 data point — open circle; asymptote from modelling 
of IRAS data ( |Serjeant & Harrison|2005| > — dotted line. Models 
are discussed in the text. Dashed line indicates the flux at which 
the integrated number density is (40 beams) -1 . 
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Fig. 2. The integrated background light at 250, 350, 500 //m from 
the HerMES counts determined in Figure [T] Dotted lines are the 
flux at which the integrated density is (40 beams) 1 . The hatched 
regions are measurements of the COBE background fLagache] 
|et al |1999[ > 



paper will be released through the Herschel Database in Marseille HeDaM 
(hedam.oamp.fr/HerMES). 
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